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ABSTRACT 

We apply a non-linear statistical method in turbulence to the cosmological perturbation theory and 
derive a closed set of evolution equations for matter power spectra. The resultant closure equations 
consistently recover the one- loop results of standard perturbation theory and beyond that, it is still 
capable of treating the non-linear evolution of matter power spectra. We find the exact integral expres- 
sions for the solutions of closu re equations. These analytic expressi ons coincide with the renormalized 
one- loop results presented bv lCrocce fc Scoccimarrol ()2006aL 120071 ) . apart from the vertex renormal- 
ization. By constructing the non-linear propagator, we analytically evaluate the non-linear matter 
power spectra based on the first-order Born approximation of the integral expressions and compare it 
with those of the renormalized perturbation theory. 

Subject headings: cosmology: theory-dark matter-large-scale structure of universe 



1. INTRODUCTION 

Cosmology now enters the era of precision cosmology. With large data set from the precision measurements of the 
cosmic microwave background anisotropics as well as the m atter density fluctuations in the large-sca le structure, the 
standard cosmological model has been fully established fe.g.. lSpergel et al]|2006l : [Tegmark et al.|[2b06[ ). The associated 
cosmological parameters are well-determined with an accuracy of less than 10% level. With improved sensitivity and 
higher precision of future observations, the modern picture of the universe will be further reinforced and one can even 
explore a tiny signature of new physics beyond the standard cosmological model. 

In fact, several ambitiou s missions for galaxy redshift survey are planed in order to reveal the nature of dark energy 
(e.g.. lAlbrecht et al.ll200^ iPeacock et al.l 12006 ^. and references therein). Among these, Wide-field Fiber- fed Multi- 
Object Spectrograph (WFMOS) may be one of the best facility capable of achi eving the percent-level mea. surement of 
baryon acoustic oscillations (BAOs) imprinted in the matter power spectrum (|Meiksin. White fc Peacoc3 [l999). The 
recent observations from Sloan Digital Sky Survey (SDSS) and Two-degree-Field Galaxy Redshift Survey (2dF GRS) 
showed that the characteristic scale of B AOs can be used as the cosmic standard ruler to determine the distance - 
redshift relation of high-redshift galaxies ()Cole et al.ll2005t lEisenstein et al.ll2005t iHiitsil 120061 : iPercival et all 120071 ). 
Since the distance-redshift relation is sensitive to the co smic expansion history, details of the accelerated expansion 
can be clarified from the accurate measurement of BAOs (jSeo fc Eisensteinll2003l ). With percent-level measurement of 
characteristic scale of BAOs, the determination of dark energy equation-of-state, parametrized hy w = Pde/pde, will 
achieve a few percent accuracy, where Pde and pde are the pressure and the energy density of dark energy, respectively. 

On the other hand, pursuit of the nature of dark energy highlights various fundamental problems which are po- 
tentially crucial for the accurate determination of dark energy equation-of-state. For example, the observation 
of BAOs requires high-precision theoretical template for the matter power spectrum in the relevant wavenumber, 
k ~ 0.1 — 0.3/iMpc~^. To achieve the required accuracy in the determination of w, several systematic effects must be 
incorporated into the theoretical predictions. Among known systematic effects, the non-linear gravitational clustering 
is one of the most fundamental building block in the theory of structure formation. The recent A^-body simulations 
' showed that the non-linear growth of matter distribution significantly alters the shape of the power spectrum and 
the acoustic sign ature of BAOs tends to be erased, where linear theory prediction of matter power spectrum is no 
longer vahd (e.g.. iMeiksin. White fc Peacockl[l999l : ISeo fc: Eisensteinll2005l ). To tackle the issue, the pe rturbation the- 
ory for gravitational clustering has been revived and has been applied to the study of BAOs (e.g.. ISuto fc Sasakil 
1991'; 'Makino, Sasaki &_Sutdll993 iJain fc Bertschingerj[l99l IScoccimarro fc Friemanlll996t iJeong fc Komatsull2006f 
Nishimichi et al. 2007). The inclusion of leading-order correctio n to the non- linear cluste ring effect somehow improves 
a performance and reproduces the A^-body results very well ijJeong fc Komatsul l2006f ) . At lower redshifts z < 2, 
however, next-to-leading order effect becomes important and the theoretical prediction with leading-order correction 
is insufficient to reproduce the A^-body simulations. 

Going beyond the perturbation theory, existing theoretical tools dealing with the non- linear gravitational clustering 
are th e A-body simulation and the fitting formula for matter power spectrum (e.g^Pcacock & Dodds 1996; S mith et al.l 
|2003| ). as well as the phenomenological approach based on the halo model (fCoorav fc ShethI I2002I . for a review). 
Currently, however, none of the reliable methods to ensure the percent-level precision exist. While the A^-body 
simulation has a potential to provide a high-precision prediction, at present, one cannot blindly trust the A'^-body 
results unless reliable and comparable counterpart is established and is fully reconciled with A^-body results. In this 
respect, development of new analytical method beyond the perturbation theory is necessary and essential for the 
progress on precision cosmology. 

In this paper, we present a non-linear statistical method to predict the time evolution of matter power s pectrum. 
Very recently, there appear several works on the statistical treatment going beyond the perturbation theory (jValageasI 
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20071 ). Based on the field-theoretical approach, the perturbation theory has been reformulated by improving the 
su nimation of the naive p erturbative expansion. The so-called renormalized perturbation theory (RPT) developed 
by ICrocce fc Scoccimarrol (|2006a) seems viable theoretical tool alternative to the A^-body simulation, suited to a 
high-precision prediction. Us ing RPT, some attempts to predict the non-linear evolution of BAOs has been reported 
(jCrocce fc Scoccimarrol l2007j ). Here, we consider an alternative statistical method accepted widely i n the subject of 
statistical theory of turbulence (e.g., lKraichnanlll959l : iLeslidlTOTSt lKanedalll981l : iKida fc Gotolll997[) . In contrast to 
the sophisticated treatment based on the field-theoretical approach, our method is rather primitive in the sense that 
the effect of non-linearity of matter power spectrum is simply described by a systematic expansion and a truncation 
of the naive perturbation. After applying the so-called reversed expansion, the perturbative expansion is effectively 
re-organized and a class of higher-order corrections is systematically re-summed (.WYld.l961i) . With this treatment, 
some non-perturbative effects are also incorporated. We derive a closed set of moment equations characterizing the 
non-linear evolution of power spectra. The resultant evolution equations consistently recover the leading-order re- 
sults of standard perturbation theory (one- loop PT). The solutions for these closure equations have exact integral 
expressions, which coincides with the one-loop results of RPT. Constructing the non-linear propagator, we attempt to 
evaluate the non-linear matter power spectra analytically. Based on the first-order Born approximation of the exact 
integral expressions, we find that the power spectra from the closure theory reasonably agree with those from the 
RPT. 

The outline of the paper is as follows. In Section|2l governing equations for the dynamics of cosmological gravitational 
clustering is presented on the basis of the fluid description of the Vlasov equation. The closure problem, i.e., theoretical 
issue on the self-consistent treatment of the hierarchy of moment equations, is briefly mentioned. In Section [31 non- 
linear statistical method in turbulence called direct interaction approximation is introduced and is applied to the 
present cosmological situation. Then, a closed set of evolution equations for matter power spectrum is obtained by 
the systematic perturbative expansion and the so-called reversion procedure. In Section 51 important properties of the 
resultant closure equations, i.e., recovery of standard PT and exact integral solutions, are discussed, based on detailed 
mathematical calculations presented in Appendices lAl and IbI In Section |5l analytical treatment for our closure system 
is presented for illustrative purpose. Constructing approximate solutions for non-linear propagator, power spectrum of 
density fluctuations is calculated based on the Born approximation of the integral solutions. The results are compared 
with those obtained from RPT, particularly focusing on the non-linear evolution of BAOs. Finally, Section|6]is devoted 
to discussion and conclusions. 

2. PRELIMINARIES 
2.1. Basic equations 

Throughout the paper, we consider the evolution of mass distribution in the flat universe, neglecting the tiny 
contribution from the massive neutrinos. To evaluate the non-linear growth of the density perturbations, hydrodynamic 
description of mass distribution is useful. Strictly speaking, this treatment is not exact and is often called the 
single-stream approximation of the Vlasov equation. However, at least in a statistical sense, it would be the best 
approximation if the scale of our interest is sufficiently large so that one can safely ignore the effect of shell-crossing. 
Denoting the fluctuation of mass distribution consisting of the cold dark matter and the baryon fluid by 5, we have 

f + -v{(i + 5M-o, (1) 

at a 

dv I 1 

— +Hv + -{v-V)v^ — V(/), (2) 
at a a 

^V^^^AnGp^S (3) 

Here, a is the scale factor of the universe, H is the Hubble parameter and pm is the ho mogeneous ma ss density field. 
Assuming the irrotationality of the fluid flow, the above equations can be recast as (e.g.. lTaruvall2000l ): 

^+He + -ViS-v) = 0, (4) 
at a 

30 1 13 

— + -{1^3wn^)H9+ — V-{{v- V)v} + -H{1 - n^)S = 0, (5) 
at 2 a'^H I 

where the quantity is the velocity divergence defined by 

V • V 

The quantity 17^, is the density parameter of dark energy satisfying the equation of sate, Pdc = wp^c, defined by 
r2^, = 8tt G pdc/{3H^). Note that the relation -1- il^ = 1 holds in the flat cosmology. 

To treat the non-linear evolution of the matter power spectrum, the Fourier representation of equations ^ and ([5|) 
is useful. To do this, we introduce the Fourier transform of the perturbed quantities: 

'5(^;^)- / 7#^e-'^-^5(fe;t), 9{x;t) ^ [ e"' ^(fe; t). (7) 
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Closure theory for cosmological power spectra 
The assumption of irrotational flow implies 



Then, the fluid equations ^ and ^ can be written as (e.g.. lTaruvall2000f) : 



H-'^^^+e{k;t)= - J ^a{k\k-k')5{k-k'-t)e{k'-t), (9) 



^^,mkpi ^ 1^^ _ 3wn^)eik;t) + ^(1 - n^)S{k;t) 



1 r fi^k' 

' (5{k' ,k-k')e{k~k'-t)e{k'-t), (10) 
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where the kernels in the Fourier integrals, a and (3, are respectively given by 

^/ ; ^ 1 , '^I'^s an 1 \ (^1-^2)1^1 + ^21^ 
a(fci,fc2) = l + ^, /3(fci,fc,) = ^-^^^-p . (11) 

For later analysis, it is convenient to introduce the vector-field notation and rewrite the equations ^ and (jlOp in 
more compact form. For this purpose, we define 

^<^^^-^')-{-etSi)/m)^ (12) 

where the function fit) is given by f{t) = d\nD{t)/d\na and the quantity D{t) being the linear growth rate. Then, 
in terms of the new time var i able rj = In-D(f) , the evolution equation for the vector quantity $a(fc;t) becomes (e.g., 
ICrocce fc Scoccimarroll2006al: IValagea£ll2007af ): 



d 

Sab ^ + ^abiv) 



/(Pk (Pk 
(27r)3 ^^^^ labc{ki,k2) Mki;t) <^c{k2;t), (13) 

where i5_d is the Dirac delta function. Here and in what follows, we use the summation convention that the repetition 
of the same subscripts indicates the sum over the whole vector components. The time-dependent matrix flabiv) is 
given by 

^abiv) = ( _ J_n - n ) _ 1 V (14) 

Each component of the vertex function jabc becomes 

ra(fc2,fci)/2;(a,5,c) = (1,1,2) 
^ (L. r. N J a(fei,fc2)/2;(a,6,c) = (1,2,1) , . 

7a6c(fci, fe2j - < p(^k,,k2)/2 ; (a, 6, c) = (2, 2, 2) ' ^^^> 

[ ; otherwise 

Note that the vertex function jabc has the following symmetric properties: "tabc{ki,k2) = 7006(^2,^1), 
7Q6c(-fei, -^2) = 7Q6c(fci,fc2), 7a6c(fei,-fc2) = labc{-ki,k2) and '-fabc{k,-k) = 0. Equation dTS]) with IE]) 
is the basic equation for our subsequent analysis. 

2.2. Moment equations 

Before addressing the non-linear statistical method, it is worthwhile to mention the closure problem for the dynamics 
of statistical quantities. 
First of all, we define the two kinds of the power spectra for fluid field Ui{k; t): 

($,(fe; rj) $b(fc'; ry)) - {27rf Soik + k') Pab{\k\; ??), 

{^a{k;ri)Mk';v')) = {'^T^fSD{k + k')Rab{\k\-T^,ri'), {v > v'), (16) 

where the bracket (■) stands for the ensemble average. The quantity Pab is the ordinary power spectra which we are 
interested in and we have the symmetry, Pab = Pba- On the other hand, the quantity Rab represents the cross power 
spectrum between different times, which will be later important when we derive a closed set of equations. Note that 
Rab ^ Rba, in general. 

Since we are specifically concerned with the time evolution of statistical quantities Pab{k] rj) and Rab{k] rj, rj'), rather 
than focusing on equation (|13p . it seems convenient to treat the moment equations for these quantities. To derive the 
moment equation for Pab, we first note that 

I ,)) ^ ,)) + ,) ^5^) , (17) 
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With a help of equation (fT3| . we eliminate the time derivative d^a/drj. Then, we obtain 

^abcd 

(5i3(fc' - fci - fe2) 7bpg(fci, fe2)^*Q(fe; v) *p(fei; f?) *g(fc2; »/)^ 

-fo(fc-fci -fc2)7ap9(fei,A;2)/-I>6(fc';?7)$p(A;i;?7)$g(fc2;?7)) ■ (18) 



Similarly, moment equations for Rab becomes 

^^^^3 - fci - fc2) 7ap,(fci, fc2)($c(fc'; r,') a>p(fci; ry) $,(fc2; 77)). (19) 

Here, we have introduced the two kinds of operators, Sabcd and h.ac'- 

^ d ^ d 

'^abcdiv) ^ ^acSbd + Sac^bdiv) + Sbd^aciv), Aabiv) = ^ab + ^abiv) ■ (20) 

Oil Oil 

Equations (fT8|) and (|19|) are not yet closed because they contain the higher-order correlation functions (or bi-spectra) . 
In order to obtain the closed set of evolution equations, it is necessary to derive the evolution equations for higher-order 
correlation functions. However, the repetition of this treatment produces an infinite number of evolution equations and 
one cannot obtain a closed set of equations. This is the so-called closure problem for dynamics of statistical quantities. 
Note that the closure problem considered here is very close to the concept of BBGKY hierarchy, but slightly different 
in some sense. The BBGKY hierarchy arises from the many-body system characterized by the Loiuville equations 
and it also appears in the linear system. On the other hand, the origin of the closure problem essentially comes from 
the non-linearity of equation (|13p . Hence, to derive a closed set of moment equations, one must devise to introduce 
some truncation procedures by approximately treating the non-linear interaction in a self-consistent manner. The 
self-consistent truncation procedure is referred to as the closure theory (or closure approximation) in the statistical 
theory of turbulence and various closure theories have been so far exploited. In what follows, we especially consider 
the so-called direct-interaction approximation as one of the reliable closure theories. 

Before closing this subsection, we define the propagator Gab{k,r]\k' ,r]'), which will later play a key role in deriving 
a closed set of equations: 

G.,{k,rj\k\,')^§§4- (21) 
d^b{k ;r]') 

where S stands for a functional derivative. It represents the influence on $a(fc;?7) at time i] due to an infinitesimal 
disturbance for $fc(fc ;??') (77 > ry'). Taking a functional derivative of equation p3p . we obtain the governing equation 
for the propagator as 

Aabiv)Gbc{k,r]\k',r]') ^2 J — ^^-^ — (5r,(fc - fei - fc2) 7ap<;(fei, ^2) 

X %iki;r^)G,c{k2,v\k',v') (22) 

with the boundary condition: 

Gabik, ri'\k\ 77') - 6ab Soik - k'). (23) 

3. CLOSURE THEORY 
3.1. Direct-interaction approximation 

In the subject of fluid mechanics, statistical characterization of turbulence for incompressible fluid flows is one of the 
major goals to understand the non-linear dynamics of Navier-Stokes equations. Among various attempts o const ruct a 
statistical theory o f turbulence, there are approaches on systematic renormalized expansion (e.g.. lKraichn an 1959; Wvldl 
119611: lLesliellT973[ ). Direct-in teraction approx i mation (DIA) is one of the best-k nown approximations and provides a 
simple truncation procedure ()Kraichnanlll977l : [Kanedairi981l : lKida fc Gotolll997f) . DIA has several desirable properties 
such as the local energy conservation and the realizability of energy spectrum. Also, the agreement with numerical 
simulations of isotropic turbulence is excellent even at high Reynolds number. Although we are especially concerned 
with the dynamics of compressible and irrotational fluid flow, the non-linearity in the cosmic fluid system essentially 
come from the same advection terms as in the Navier-Stokes equations. In this respect, DIA is a promising method 
to give a quantitative prediction of matter power spectrum. 

To derive a closed set of evolution equations in the DIA, one transparent and intuitive way is to decompose the true 
field into the direct-interaction (DI) and the non-direct-interaction (NDI) parts. Let us consider particular Fourier 
modes (fc, p, q) and especially focus on the time evolution of '^a for a specific Fourier mode fc. Through the non-linear 
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term in the right-hand side of equation (|13p . the time evolution of $a(fc; v) is determined by the infinite sum of three 
Fourier modes. We can formally decompose the quantity ^a{k; rf) into 

r?) = ^ ^(Di)(^. (24) 

In the above expression, the quantity $1^^^^ denotes the DI field, whose time evolution is determined by the direct 

interaction with the particular Fourier modes, p and q. On the other hand, the quantity is defined as a fictitious 

field without the direct interaction between three modes A;, p and q. 
In the system governed by dimensionless equation (fT3|) . the strength of the non-linearity is characterized by the 



number of interacting Fourier modes. In the non-linear regime, one naturally expects that NDI field ^i^'^^' plays the 
most dominant part in the time evolution of the quantity ^a{k] rj). Hence, in the DIA, we treat the DI field as a small 

perturbed quantity, relying on the assumption ^'i'^^'^ » In addition, we put the following assumptions: (i) 

Gaussianity of the NDI field; (ii) statistical independence among the modes without direct interaction; (iii) statistical 
independence between the NDI field and the propagator. These assumptions basically come from the physical intuition 
for fully developed turbulence. In the presence of the infinite sum of the quadratic interaction, the fluid fields are 
expected to be nearly Gaussian, as naively indicated from the central-limit theorem. Also, it seems plausible that the 
initial correlation between different modes or quantities tends to be lost and the fluid quantities become statistically 
independent along the course of the non- linear interaction. Then, relying on these assumptions, we systematically 
expand the moment equations and the governing equation for propagator. Evaluating the ensemble average by using 
the for mal solution of DI field in terrns of the NDI fields, the closed set of equations for NDI fields can be finally 
derived CKraic hnanl[l959l : IKida &: Goto! [l997l : [Goto fc KidallfOOl . 

In t he following, we shall derive a closed set of equations in an alternative route. As it has been shown bv lGoto fc Kidal 
(1998), the resultant closure equations by DIA is identical with those o btained fro m the so-called reversed expansion 
procedure by introducing the fictitious parameter^ (see also lLeslidll973t lKraichn;aiilll977t lKanedal[T98I ). The reversed 
expansion procedure seems rather straightforward and the assumptions made in this procedure are relevant for the 
present cosmological situations. 

3.2. Derivation 

Let us first introduce the fictitious parameter A, which represents the strength of the non-linearity, and consider the 
weakly non-linear regime. We rewrite the field $a and the propagator Gab with 

$a(fc;77) = Aia(fc;r,), Gab{k,v\k' ^ ^Gab{k,v\k' ,fi'). (25) 
In terms of these, the basic equations (fT3|) and ((22)) respectively become 

Aafc(?y) $6(fc;?7) = A J ^^^^^ fo(fc - fci - fc2) 7afcc(fci, fc2) $fc(fci; ??) $c(fc2; ??), (26) 
^ ~ f (Pki (Pk2 

Aafc(?y) Gbc{k,'n\k\T]') ^2 X J — (27r)3 ^D{k - ki - k2)japq{ki,k2) 

X $p(fei;?7)G,,(fe2,77|fc',?7')- (27) 

Our task is to derive the consistent closure equations from the moment equations (fT5)l and ([TO]) with a help of the 
factitious parameter A. To do this, we regard A as a small expansion parameter (i.e., A ^ I) and put the following 
assumptions: 
Assumptions 

(i) The evolution of field ^aik; r/) is started from a tiny fluctuation and the development of non-linearity is mild. 

(ii) The field ^a{k; rj) is a homogeneous random field and the statistical property of $a(fc; v) is approximately described 

by the Gaussian statistics. 

(iii) At the leading order, the field $q and the propagator Gbc are statistically independent from each other. 

Then, based on the perturbative calculation, we first evaluate the higher-order correlation terms in the moment 
equations (jlSp and (|I9p . Inverting the perturbative expansion by a formal replacement of the perturbed quantities 
and setting the fictitious parameter to A = I at a final step, we obtain a closed set of evolution equations. The set 
of equations derived here may be regarded as the result of renormalization and/or resummation of the perturbative 
expansion, which will be applicable to the nonlinear regime of gravitational clustering beyond the standard perturbation 
theory. 



^ In subject of turbulence, this fictitious parameter corresponds to the Reynolds number. 
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3.2.1. Naive perturbation 

Based on the assumption (i), let us evaluate the nonlinear terms by a perturbative expansion of the small parameter 
A. To do this, we expand the quantities and Gab as: 

Gabik, v\k', v') = G'-ab (k, vW, ?/) + A (fc, rj\k', ,,) + •••. (28) 
The boundary condition for propagator in each order is 

Gf^{k,7^'\k',v') = SabSD{k-k'), G^2{k,^'\k',v) = 0. (29) 

Also, the boundary condition for the first-order quantity ^i^'' at the initial time 770 is 

$«(fc;77o) = 0. (30) 

Since the equations for zeroth-order quantities are source-free, no mode-mode coupling occurs. Thus, the zeroth-order 
propagator g'"^^ satisfying the boundary condition ([^5]) can be expressed in the following form: 

&;^^{k,T^\k','n')^Gab{k\ii,i) So{k-k'). (31) 

Note that the function Gabikli],?]') is a dimensionless quantity. Using this functional form, we formally write down 
the first-order solutions to the quantities ^l^'' and G^^^ From ([26|l and ([27]) . we obtain 

d^kid^k2 

ari Uac[i«\ri,'n ) / 

'm 



$W(fe;r;)=y dry" G,,(fc|ry, 77") / (2^)3 Soik-k.^k,) 



X 



7cp,(fei , fe2) (fei; 77") $ W (fe2; 77"), (32) 



G«(fc,77|fe',ry') = 2£d,7" Gac(fckr;") 1 '^^^0^ 5D{k - k, - k,) 

X 7cp,(fei,fc2)$f (fci;r,")Gj)(fc2,,7"|fc',,7') (33) 

for > rj' > rjQ. 

According to the assumption (ii), one may treat ^i^"* as a Gaussian random variable. Then, applying the perturbative 
expansion (|28|). we calculate the lowest-order non- vanishing contribution to the three-point correlation in equations 
(fT8|) and ([T9|) . Let us first deal with the right-hand side of equation (fT8|) . At the lowest-order contribution, the 
three-point correlation becomes 



^$,(fc; ,7) $p(fci; 77) $,(fc2; r/)) ^ A ($i")(fc; r;) $(")(fci; 77) $^'^^2; ^7)^ 

+ A($(0)(fc;7;) $W(fci;77)$(°)(fc,;^); 

+ A ($W(fc;77)$(°)(fci;77)$W(fc2;^)), (34) 
which are of the order of C'(A^). Substituting the formal solution ([5^ into the above, we obtain 



^aik;r])^p{ki;r])<S>g{k2;v)) - A /" dry' /" 



(2^)3 

^D(fe2 - P - q) G,z(fc2|77, 77') ^irs{p, q) ^) ^) 51°' (p; ^') ^i^'C?; ^') 

+<5z5 (fei - p - q) Gpz (fci ,7') 7/.. (P, q) (5i°' (fc; 77) $(0) (p; 77') (q; 77') (fc2; 77) 

ik-p-q) Gal {k\v, v') iirs {p, q) {^'^ (p; v') $i°' (g; V) $f (fci ; ^) (fc2; r/) 



(35) 

The above expression is further reduced if we use the perturbative expression for the power spectra (jl6p : 

Pab{k;v)^ X' Pab{k;v) - {p^°\k;v) + 0{X^)] , 

lab{k;v,v')= X' Rab{k;v,v')-X' {i?i",'(A; 77, 77') + ^(A^)} . (36) 



Ra 
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The leading terms and Rf^j are defined by 

r^) ry)) = i2n f 5n{k + k') P^l\k; ry), 

;$(°)(fc;,7)$i°)(fc';ry')) - {2n f Sn{k + k') R^°,\k;rj,7^'). 
Then, after some algebra, equation ([55)1 is finally reduced to the following form: 

%aik; rj) $p(fei; ?7) $,(^2; t)) ^ (2^)^ Soik + fei + fea) A ^^i^^ (fc, fci, fca; v)- 



(37) 



The function Fipg is explicitly written in terms of the two-time correlation i?^^'' : 



F^ll{k,ki,k2;v)^2 I dii [Ggi{k2\v,v')jirs{k,k,)Rl^J{k;rj,7^')R^^Xki;V,v') 



(38) 



In deriving equation psp . we have used the symmetric properties of the vertex function, i.e., ^abc{ki, k2) — 7acfc(fc2, ^i) 
and jabc{k, -fc) = 0. 

In similar manner, we perturbatively evaluate the three-point correlation in equation (jl9p . The resultant expression 
becomes 



(39) 



$,(fc'; 77') $p(fei; 7?) $,(fc2; r,)) ^ (2^)^ Soik' + fci + ^2) A if (fc', ^1,^2; r,, t?'), 

(?7 > 77') 



with the function K, being 



bpq 



drj" {i?(0)(fc';r/, '?")©('/ - v") + Ri^Hk' ;v" ,v'Mv" ^ v')} 
G,i{k2\v, V") lirsik', fci) Rf}{ki;Ti, ri") + Gpi{ki\rj, rj") j^k' , k^) R^^J{k2;v, v") 
dfj" Gci{k'\fi', v") jirs{ki,k2) R^"}{ki;f], 77") R^°Xk2;v, v")- 



(40) 



Here, the function (d(t) is the Heaviside step function 

Now, summing up the perturbative expressions of three-point correlations, the moment equations (llSp and ([1 
respectively become 

lbpq{q,-k - q) F^l\{k,q,~k ~ q] r]) 



-'abed 



{r,)P,d{k-r,)~\' 



{2n) 



'lapq{q,k~q)Fl^l{^k,q,k-q:r^) + 0{X% 

d^q 
(2^ 

3.2.2. Non-linear propagator 

Next, we evaluate the higher-order terms in the governing equation for propagator (j27p . To do this, we first notice 
that the propagator Gab is no longer deterministic. Because of the interaction with the random field $q, the time 
evolution of Gab also exhibits stochastic nature. We thus treat the equation (l?7)) statistically: 

d^ki d^k2 



Aabiv) Rbcik; 7?, rj') j ^apqiq, k q) Kgl{-k, q,k q; rj, 77') + ©(A^ 



(41) 
(42) 



Aah(r/)(G;,e(fe,7;|fe',77'))=2A 



{2n) 



3 Soik - ki - k2)^apq{ki,k2) 

X (^p{ki;ri) Gqcik2,v\k' 



(43) 



The higher-order term in the right-hand side of equation (|43p is perturbatively evaluated as 

%{k^:,Tj) G,,(fc2, v\k', v')) - (^f iki;v) &qj {k2, v\k' , v')) 

+ \ ($«(fci;7;)G(°)(fe2,77|fc',77') 
-l-A ($f(fci;ry)GW(fc2,?7|fc',??') 



(44) 



8 Taruya & Hiramatsu 

In the above equation, the first two terms at the right-hand side become vanishing because of the assumption (iii). 
The only non- vanishing contribution comes from the last term, which can be recast as 

($f (fci;,7)GW(fe2,77|fc',r/))=2 J\r^" G,i{k2hv") J '^0^ S^ik^ - p - q) 

X 7/..(p, q) (fci; 77)$(°)(p; v")) Gf}{q, v"\k',v') 

= 2 5D{k' -k,-k2) Pdv" Ggi{k2\v,v") 

Jri' 

X7/,,(-fci,fc') &^){ki;T^,ri") Gsc{k'\Tj,ij"). (45) 

In the last equality, we have used equation ([3T|) and the definition of Rf^^ (see Eq. [36] below) . Hence, the perturbative 
evaluation of equation (|27p becomes 

^ /- , \ / n f d^q I 

A^ab{-n){Gbc{k,il\k',-q')j^ 'iX^ Soik-k') J d-q" j Japq{q,k - q)'-firs{-q,k') 

X Ggi{\k-q\\rj,fj")R'^J{q\;V,v")GUk'W\v') + O(X'). 

(46) 

up to the contribution of 0{X^). 

3.2.3. Reversed expansion 

We are in position to employ the procedure of the so-called reversion to rewrite the perturbative expressions (|4ip . 
dm) and (USD. 

Recall from the perturbative expansion ((36)) that the 0{X'^) and the higher-order terms of the power spectra can be 
expressed in terms of the P^^-* and Rf^^ in principle, since the formal solution of higher-order quantity is always 
written in terms of with a help of propagator of zeroth-order, Gf^. Also, the ensemble average of the higher-order 
propagator is formally expressed in terms of the zeroth-order quantities, G^^fj and R^^i^ ■ We then regard the expansion 
(j36p as equations for Pab and Rab, the solutions of which are written in powers of A as 

Pa°\k;v) = Pabik;v) + 0{X^), R^abik;V,v') = Rabik;v,v') + 0{X^). (47) 
Similarly, we may write 

5D{k - k') Gab{k\ii, li) - ICabik, v\k', t]')) + ^(A^). (48) 



This procedure is called the reversion and it corresponds to the resummation of the perturbation series. Thus, at the 
leading order, equations ([4T|) . (|42|) and (|46)) are written in terms of the true field variables Pab, Rab and Gab- 
Now, we do not necessarily treat A as the small parameter. This is just a book-keeping parameter and we finally set 
it to unity. Dropping the tilde over the quantities Pab, Rab and Gab, we at last reach the closure equations: 

- f d^q r 

'^abcdiv) Pcd{k;v) ^ J j^-j^ [lbpq{q,-k - q) Fapq{k,q,-k - q; ri) 

+Japq{q,k~q)Fbpg{-k,q,k-q;ri) , (49) 
Aab{v)Rbcik;r],'q') ^ J -fapq{q,k - q) Kcpq{-k,q,k - q; 7],!^'), (50) 

r f d^q 

Aab{v)Gbc{k\f],'n') = 4 J di]" J j^—^1apq{q,k-q)-firs{-q,k) 

xG,i{\k~q\\7^,7^")Rpr{q;v,v")GsM',v')- (51) 
The explicit expressions for the kernels Fapq and Kcpq are summarized as 

r r 

Fapq{k,ki,k2;r]) ^ 2 / dr]" 2Gqi{k2\r],ri")jirsik,ki)Rar{k;'q,r]")Rps{ki;r],ri") 

+ Gai{k\r],ri")-firs{ki,k2) Rpr{ki;f],T]") Rqs{k2;v,v") , (52) 

Kcpq{k',ki,k2;r],ri') 

= 4 drt" Ggi{k2\v,v")lirs{k',ki)Rp,{h;f],7^") 
•'no 
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X r?', 77")e(77' - ,7") + i?,,(fc'; ,7", r^')Qir^" - ,/)} 

" drj" Gel (k'W, v") lirs (fci , fc2) Rpr {h ; r/, ri") R^, (fcz; r/, 77"), (53) 

where we have used the fact that the functions Fapq and K^pq always appear as the product of ^bpqFapq and japqKcpq- 

4. PROPERTIES OF CLOSURE EQUATIONS 
In this section, we discuss some important properties of the closure equations derived in previous section. 

4.1. Recovery of one-loop perturbations 

The closure equations in previous section have been derived in somewhat non-trivial manner by the reversed ex- 
pansion procedure. While we employ the perturbative approach when evaluating the moment equations, the final 
governing equations for matter power spectrum become the non-linear coupled system and it seems unclear whether 
the power spectra calculated from the closure equations consistently recover the results of the standard perturbation 
theory at some level. In this respect, the recovery of the perturbation calculation may be a fast important check for 
the usefulness of the closure approximation. 

In the standard treatment of the perturbation theory, the field <I>a is assumed to be a small perturbed quantity and 
is expanded as 

*a = + + $1^) + • • • . (54) 

Substituting the above expansion into the evolution equation (fT3|) . we systematically derive the perturbation equations 
and through the order- by-order treatment, the solutions for higher-order quantities are expressed in terms of the linear- 
order quantity ^i^-* . Further assuming the Gaussianity of the linear-order quantity <I>i^^ , the power spectra can be 
summarized as 

Pabik) = P!.l'\k) + {Pif (fc) + Pif (fc)} + • • ■ , (55) 

where the first term in the right-hand-side of equation is the linear power spectra and the quantities in the curly 
bracket is the so-called one-loop corrections to the power spectra, given by 



$ii)(fc; 77)$«(fe';,7)) = {2^)^ Soik + k') P^J,'\k-^i), 
<i>i2)(fc;,7)$f = (2^)3 ^^(fc + fc/) Pif (fcjry), 

$(1) (fc; 77)$f ) (fc'; ^) + $(3) (fe; ry)$« (fc'; ,7)) = {2irf Sn {k + k') P^^'^ (fc; ry). (56) 
In Appendix |A1 we show that the linear plus one- loop power spectra satisfy the following evolution equations: 

^a,cM {Pif\k) + Pif\k) + Pif (fc)} 

{2-kY |'^^(^ + ki + k2) 76p<j(fei, ^2) Fapq{k, ki,k2; 77) 

+ 6D{k ki k2)^apq{ki,k2) Fbpq{-k,ki,k2;ri)Y (57) 

where the function Fapq exactly coincides with the definition (|52p , with the replacement of the non- linear propagator 
and power spectra. Gab and Rab, with those of the linear counterparts, gab and R^^b"^' 

^abiv) gbciv, v') ^ac Sd{v - v') (58) 

and 

($(1) (fc; ry)$« (fc'; r/)) = (2^)^ fe(fc + fc') R(^^'\k; 77, 77'), (77 > 77'). (59) 

Thus, in the weakly non-linear regime, the closure approximation faithfully recover the one-loop results of the 
standard perturbation theory. This will be manifestly apparent in next section by calculating the power spectra. A 
great emphasis is that among several non-perturbative approaches, the closure equations as non-linear coupled system 
also have the potential to go beyond the perturbation theory. In fact, our closure theory is bas ically equivalent to the 
the 2PI effective action approach by 'Valagea^ ()2007al ) and the one-loop level of the RPT by ICrocce fc Scoccimarrol 
([2006a ) . This will be explicitly shown in next subsection when we obtain the exact integral expressions. 
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4.2. Exact integral solutions 

The closure equations derived before seem rather complicated and analytically intractable because of its non-linearity 
and non-locality. In practice, a sophisticated numerical treatment is required to get the exact solutions for closure 
system. Note, however, that the closure equations possess the exact integral expressions for the power spectra Pat and 
Rab, which are formal solutions of the closure equations: 

Pab{k; ri) = Gac{k\v, Vo) Gbd{k\v, Vo) Pcd{k; rjo) 

+ drii dr]2Gac{k\r],rii)Gbd{k\T],ri2)^cd{k;r]2,r]i). (60) 

Jr]„ Jr]o 

Rab{k; 77, 77') = Gacik\ri, 770) Gbd{k\ri' , rjo) Pcd{k; m) 
rn rn' 

+ drji dri2Gac{k\T],r]i)Gbd{k\v',V2)^cd{k;r)2,Vi), (61) 

The above expressions contain the function <I>(fc; 771, 772), which represents the non-linear mode-coupling between dif- 
ferent Fourier modes, given by 

$ab(fc; 7/1, 7/2) = 2 J j^-^-farsiq,k ~ q)-fbpq{q,k - q) 

X {^pr(g; 771, 772)i?qs(|fc - ql; 771, 772) 6(7/1 -7^2) 

+i?rp(g; 772, 77i)-Rs<;(|fc - q-l; 772, 771) 6(772 -77i)|. (62) 

Note that the mode-coupling function $ possesses the following symmetry: ^ab{k; 7/1,7/2) — ^baik; 772, 771). In Appendix 
[b1 the integral expressions (|60p and (|61|1 are indeed compatible with the closure equations if the mode-coupling function 
$ is given by equation ([5^ . 

The integral expressions given above have been also derived based on th e RPT by ICrocce fc Sco ccimarrol 
(^2006a^ and/or through th e path-integral formulation by IMatarrese fc Pietronil (j2007l ) (see also geas 2007ai 
ICrocce fc Scoccimarrol l2007l ). although their derivations are quite formal. In contrast to their formal expressions, 
our integral solutions have the explicit functional dependence of the mode-coupling function ^ab on the power spectra 
Rab and the propagators Gab- In the language of the RPT, this corresponds to the renormalized expressions for the 
mode-coupling power up to the one- loop order^. In this respect, the closure equations is a non-perturbative description 
of the power spectra going beyond the perturbation theory, and have an ability to predict the matter power spectra 
accurately, the result of which will be comparable to the one-loop results from RPT or the path-integral approach. 
Note, however, that our closure system has time evolution of the non-linear propagator, whose governing equation has 
been also derived by the self-consistent truncation of the higher-order corrections. On the other hand, in the RPT, no 
such truncation is considered in deriving the integral expressions. This may cause a major difference in the prediction 
of matter power spectra, which we will address in detail in next section. 

5. ANALYTICAL TREATMENT OF NON-LINEAR POWER SPECTRUM 

In this section, we compute the power spectrum of mass density fluctuation from the closure equations. Based on the 
exact integral expressions, we employ the Born approximation to obtain the analytic expressions for power spectrum. 
Further, the approximate expressions for the non-linear propagator is obtained by matching the two asymptotic 
behaviors. Combining these results, we evaluate the power spectrum of mass fluctuations and compare it with the 
one obtained from the RPT ijCrocce fc Scoccimarro 2007). In what follows, the following cosmological parameters are 
adopted to compute the power spectra : fim.o — 0.27, fJb.o = 0.043, il„,,o = 0.73, w — —1, h = 0.7, as = 0.8, and 
n. = 1. 



5.1. Born approximation 

In practice, numerical treatment to directly solve the closure equations would be essential for an accurate evaluation 
of non-linear power spectrum. Nevertheless, analytical evaluation of the power spectrum is instr uctive and very helpful 
to understand the behavior of the non-linear corrections incorporated into the closure system. ICrocce fc Scoccimarrol 
((20071) recently applied the first-order Born approximation to the exact integral expressions for the power spectra and 
derived the analytic expressions correctly up to the two-loop order. In the following, adopting the same approximation, 
we will analytically evaluate the power spectrum. 

Here, the term, Born approximation, means the iterative approximation scheme to evaluate the integral equations 
like (|60p and (|6H) . Applying the Born approximation to the integral solutions, the power spectra Pab and Rab are first 
evaluated by substituting the linear-order quantities into the right-hand side of the expressions ([60]) and ((6T|) . This 
is the first-order Born approximation and it can be further improved by repeating the iterative substitution of the 
leading-order solutions into the right-hand side of the integral solutions. This treatment can be accurate in principle 



^ Correctly speaking, closure approximation or DIA drops all the corrections arising from the vertex renormalization (jWvldl II 9611 V 
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and the higher-order corrections become neghgible as long as the contribution from the mode-couphng function is 
small, compared to the first term in right-hand side of integral equations. 

In next subsection, we will present the approximate expression for the non-linear propagator, G^j'"^™. Here, assuming 
the analytic form of Gab, we derive analytical expressions for the quantity Pab{k',ri) by substituting the iterative 
solutions of the different-time power spectra Rab{k] rj, rj'). Let us denote the linear power spectra given at initial time 
by 770). For a sufRciently small value of t]q, the late-time evolution is dominated by the growing-mode solution. 

We thus put 

P^i"(fc;%) = e2''°Fo(fc)(l }) (63) 

with the quantity Po{k) being the linearly extrapolated spectrum at the present time. Then, from equation (j6ip . 
different-time power spectra are iteratively evaluated as follows: 

Rab{k- f^, rj') = i?«(fc; 77, rj') + r';^^ (fc; ry, ,7') + • • ■ ; 
R^'hk;V,v')^Ga{k\r,,m)Gb{k\v\vo)e''"'Po{k;vo), 

R^I'{k;i,ri')^2 j ^Ia{k, q;v, Vo)Ib{k, q;v', Vo) 6^"^° Po{q) Po{\k - q\). 

with Ga = Gal + Ga2- Substituting the iterative solutions of Rab into the integral expression (|60p . one obtains 

Pabik; 77) = P^Sik; V) + Ptf\k; 77) + P^f'\k; 77) -f • • • ; (64) 



P^]!ik;v)^Ga{k\7j,rio)Gb{k\v,Tjo)e'^°Po{k), 

/73 
^ /a(fc, q; 7;, 770) Ibik, q, 77, 7;o) e^^^ Pa{q) P^Qk - q\), 

Pab^\k;Tj)=8 J J ^^Jaik,p,q;T], Tjo) Mk, p, 5:77,770) 

X p„{\k-p\)Po{q)Po{\p-q\). 
Here, the source functions la and Ja are respectively given by 

/a(fe, q;77, 770) = / dr]' GaiiM-q^v') lirs{q,k - q) Gr{q\ri' ,r]o) Gs{\k - qWri' ,r]o), 



V 



Ja(A;, p, q;77, 770) = / drji / dr]2Gaiik\r],rii)^irsip,k - p)Grcip\m,m) 

Jrin J Tin 



' Vo -J Va 

X lcpq{q,P - q) Gp{q\T]2,m)Gqi\p - q\\m,Va)Gsi\k - p\\r]i,r]o). 



Compared to the results given bv lCrocce fc Scoccimarrd ()2007t) . we find that the first and second terms in equation 
(|64|) . P^^^ and Pab \ exactly coincide with their expressions, G^Pq and P^c"^ in RPT, respectively. On the other 
hand, the third-order term Pj^]^^^ is very similar to the term P^c"'' of ICrocce &: Scoccimarrd (|2007[ ). but the factor 
2 is different. This discrepancy may be caused by the fact that the term P^'^') given above is associated with the 
higher-order contributions to the Born approximation of the mode-coupling function (|62p . while the expression for 

P^'^^ has been properly derived as the leading-order contribution to the mode-coupling term at two- loop order. For 
the influence of the higher-order corrections, we leave the discussions in section [5.31 



5.2. Approximate solution for non-linear propagator 

Having obtained the analytic expressions for power spectra, our remaining task is to get the approximate solution 
for the non-linear propagator, G^ ^^*^™, from the closure equation ((?T|). Here, we just follow the procedure suggested 
by ICrocce &: Scoccimarrd (|2006bf ) and construct the approximate solutions restricting their validity to the low-fc or 
the high-fc regions. Matching these asymptotic solutions appropriately at an intermediate regime, we obtain the global 
solutions for the non-linear propagator. 

5.2.1. Solutions at one-loop order 

Let us first consider the \ow-k limit of the non-linear propagator, where the perturbative treatment is safely applied. 
As it has been shown in section B?T1 our closure system consistently reproduces the one- loop results of the perturbation 
theory for power spectra Pab- Indeed, this is also true for the propagator Gab- Here, we explicitly write down the 
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perturbative solutions at the one-loop level, 5G^j,'°°''(fc|7], 77'). Using the linear propagator gab, perturbative solution 
of equation ([5T|) is generally expressed as 

^Gll^°°^ik\Tj,r]')=4: J dr]i gac{v,Vi) J dm J lcpq{q,k - q) 



x7/rs(-q, k) .9,; (771, 772) R^pr{q; m,m) 9sb{v2,v'), (65) 



with the quantity being the linear-order solution for different-time power spectra. In the approximation that the 
time-dependent matrix Qabiv) given by equation p4p is replaced with the constant matrix with fi^^, = and / = 1, 
the analytical solution for linear propagator q„ b satisfying the evolution equation (|58p is obtained and is given by (e.g., 
iCrocce fc Scoccimarroll2006al: IValageaj|2007ah : 

gab(Vi,V2) = [ — ^ — 132/ 5 1-3 3 j/Qw^^a) (66) 

with Q[x) being the Heaviside step function. For the different-time spectrum R^b, we have 

i?!ifc"(fc;^i,^2) = e''^+''^ (^ll^Poik), (67) 

where we have neglected the contribution from the decaying mode. 

Substituting the quantities and ([57]) into the next-to- leading order solution ([55)) . we first perform the time 
integrals over 771 and 772. As for the three-dimensional Fourier integral, we can write d^q as 2'iTdqq^dx, where we 
have performed the integral over azimuthal angle. The variable x is the cosine of the angle between k and q, i.e., 
X = k ■ q/(kq), and the integr al over x is performed analytically. A straightforward but lengthy calculation leads to 
(jCrocce fc Scoccimarro|[2006bD : 



where the matrices Xab and Yab are given by 



_ 2,' /"(^ ~ ^')-^(^) - (^aiv - ri')Kk) a{r, - r^')f{k) - (3,{r^ - T^')h{k) \ 
^ ya(77 - 77')g(fc)+7g(77 - 77')/i(fc) «(?? - ??').9(fc) - ^ 7g(?7 ^ ?7')«(^)y ' ^^^^ 

" ^ " U(77 - 77')5(fc) + - rly{k) 5(77 - 77')/(fc) - \ ^ v)h{k) ) ' 



The above expressions contain time-dependent functions a, Pgd, 7g,d and i5 and scale-dependent functions /, g, h 
and i, whose explicit expressions are summarized in Appendix[C] Note that the large-A: limit of the above functions 
satisfies 

/, 9, h, I ^-i(fcav)2 (71) 
with the quantity (Tv being the velocity dispersion for linear fluctuation defined by 



(27r)3 g2 



(72) 



5.2.2. Solutions in the high-k hmit 

Turn next to consider the high-fc limit of the non- linear propagator. In the evolution equation (|5ip. we take the 
limit k — > +00, while keeping q finite. In this limit, the vertex functions behave like 

Icpqiq, A; - g) ~ i ^-^ Scq ^2p, llrsi-q, k) ^ - ^ ^^-^ Sis 52v 

To estimate the leading-order behavior analytically, the different-time spectrum i?pr(9; J/i, ??2) in equation (jSip is also 
treated approximately by replacing it with the linear-order quantity i?p", given by equation (j67p . 
Then, the governing equation for non-linear propagator (|5ip is greatly simplified and we obtain 

Kb{l) Gbc{k\v, v') = -{kcT. f r dv" Gab{k\v, v") Gbc{k\v" , v') 6"+"", (73) 



where the quantity Cv is the rms fluctuation of the linear velocity given by equation (|72[) . To solve the above equation, 
we adopt the following ansatz: 

Gab{k\r], 77') = gabiv, ??') f{k\r] - 77') 
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Fig. 1. — Approximate solutions for non-linear propagators Gi{k\z, Zinit) = Gn + G12 {left) and G2{k\z, Zinit) = G21 + G22 (Tight) 
as function of wave number. In each panel, solid lines show the results from closure theory, while the dashed lines are the propagators 
based on the renormalized perturbation theory (RPT) (ICrocce fc Scoccimarro 2006li). For comparison, we also plot the results from the 
perturbation theory (l-loop). From left to right, the lines indicate z = 0, 2 and 5 with initial redshift zinit = 35. 

with the initial condition, /(fc|0) = 1. Using the basic property of the linear propagator, gabiVj'>l')9bciv' jV") — 
ffafc(?7,'7"): equation ([73]) is rewritten with 



d_ 
97 



fik\r)^-{ka,f / dr' f{k\T-r')fik\T'). 
Jo 



(74) 



— e'' . The above equation has the analytical 

+00 




Here, for convenience, we introduced the new time variable r 

solution (| Valageasll2007af l . Writing the Laplace transform of the function / as f{k\s) = ds e~^'^ f {k\T) , we have 

sf{k\s)-l = -{ka^f[f{k\s)f. 
The solution of this equation satisfying the limit, / ^ for s +00, becomes 

Jim 



1 



-s+ y^s'^ +4{ka^y 



2{ka^y 

The inverse Laplace transform of the above expression is well-known and can be read off from the mathematical table: 

Ji{2x) 



fik\T) 



X = ka^T 



(75) 



with the function Ji{x) being a Bessel function of the first kind. Hence, the non- linear propagator in the high-fc limit 
finally becomes 



^ ,^ f „ Ji{2ka^{e^ - )) 

Gab[k v, 77 ) = 9ab{v, V — ; — rr; itt — ■ 

kav{e^ — C ) 



(76) 



5.2.3. Matching the two solutions 

The asymptotic behaviors obtained in section 15.2.11 and 15.2.21 have an overlapping region in which both of the 
approximations are applied. Therefore, matching these two solutions, one can obtain a global solution which would be 
a good approximation for the full propagator G. Let us recall from equation (|7T|) that the propagator including the 
one-loop correction has the following asymptotic form: 



g,,(77,,/) + <5G^;'°°P(A:|r7,,y') 



where we have only considered the dominant terms at 77 — > 00. On the other hand, the propagator in the high-fc limit, 
([TG]) . is perturbatively expanded as 



I X 

Gab{k\ri,ri') ~ gab(ri,ri') {l - — + 



A:av(e''-r?''). 



Comparing these two expressions, the approximate solution smoothly matching these asymptotic behaviors at 7] ^ 00 
may be 



3Pii2Pi2 
3 P21 2 P22 



2Qii ~2Qi2 
-3Q21 3Q22 



(77) 
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Here, the matrices Pat and Qab are defined as 

^ab ~ , Qab— i'oj 

Xab Yab 

with Xab = l^Xab]^^^ and Yab = \2Yab^^'^ (see Eqs.[B2] and [7D] for definitions of Xah and Yab)- In the weakly non- 
linear regime, the propagator G^^^^°^ correctly reproduces the one-loop results. In the large-fc limit, the function (j77p 
asymptotically approaches the solution ([76]) . 

Note that the approximate propagator (|77p is derived in the same way as done in RPT of iCrocce fc Scoccimarrol 
()2006hf ). although the functional dependence is somewhat different because of the different high-fc behavior. In RPT, 
the matrices Pab and Qab defined above should be replaced with 

Pab -> exp(Xafc), Qab ^ exp(yQ&), 

which lead to the asymptotic behavior, c^'^pp'"'') _^ gafc exp(— a;^/2), in the high-fc limit. 

Figure [1] shows the propagators Gi — Gii+G'12 {le-fi) and G2 — G21+G22 {right) multiplying the factor £'(2;int)/£'(z), 
for specific redshifts z = 0, 2 and 5, and with initial redshift Zinit = 35. The solid lines represent the approximate 
solutions obtained by matching the two asymptotic solutions. While the dotted lines show th e results from the one- 
loop p erturbation, the dashed lines indicate the non-linear propagators obtained from the RPT (|Crocce fc Scoccimarrol 
l2006bD . As increasing the wave number fc, all the results exhibit the decaying behavior and the characteristic scale 
of the decay is shifted to low-fc as decreasing the redshift. A closer look at small scale (high-fc) reveals that the one- 
loop propagators, gi^2 + 2°°^; show unphysical behavior, which eventually become negative and tend to diverge. 

On the other hand, the approximate solutions G^'^f'"'^ obtained from the closure theory and RPT asymptotically 
approach zero as fc ^ 00. These damping behaviors are regarded as the non-perturbative effect as a result of the 
renormalization and/or self-consistent closure, which effectively takes account of the infinite series of higher-order 
corrections. Nevertheless, there exist some differences in the damping behaviors of the propagators. While the 
propagators in the RPT exhibit the exponential damping, the approximate solutions in closure theory show a damping 
oscillation. These differences may affect the final result of power spectrum. This point will be carefully discussed in 
next subsection. 

5.3. Results and discussion 

Having provided the basic ingredients for calculating the power spectrum, we now present the analytic results for 
power spectrum of density fluctuations, i.e., Pii(fc;z), and compare those with the results obtained from the RPT, 
particularly focusing on the characteristic scale of the BAOs. 

Figure [2] illustrates the overall behaviors of non-linear power spectrum of density fluctuations P{k;z) = Pii(fc;z) 
given at z = 1, based on the Born approximation (|64ll . Here, the contributions to the total power spectrum up to the 
first-order Born approximation, i.e., P^^^(fc) and P^^^^(fc), are separately plotted. Thin and thick lines represent the 
res ults from RPT and the closure theory, respectively. The result from RPT is basically the same one as presented 
bv ICrocce fc Scoccimarrol ([2007), although they further considered the higher-order contribution coming from the 
two-loop correction. Due to the damping behavior in the non-linear propagators, each contribution to the total power 
spectrum rapidly falls off in both predictions and their amplitudes become significantly lower than the linear power 
spectrum on small scales (labeled by Piinear(fc)), where the differences between the two predictions become manifest. 
Turning to focus on the scales larger than the damping scale, the contribution coming from p('^) becomes maximum 
around fc ~ 0.2h Mpc~^, where the non-linear enhancement of the power spectrum P{k) can be seen and the differences 
between closure theory and RPT become fairly small. In particular, for the scale of our interest on BAOs around 
fc ^ 0.01 — 0.3ft. Mpc~^, one cannot clearly distinguish between both predictions from Figure [21 

To enlarge the differences between the two predictions and to clarify the non-linear behaviors of the BAOs, in Figure 
[3l snapshots of the power spectra divided by the smooth linear spectrum, P(fc) /Pno-wiggie(fc) are plotted for specific 
redshifts z = 3, 2, 1 and 0.5 together with the A^-body results kindly provided bv lJeong fc Komatsul ()2006l ) (except for 
z=0.5), while in Figured! we present the logarithmic derivative of the power spectra, d In P(fc)/(i In fc. All the results 
are plotted in linear scales. The smooth power spectra, Pno-wiggic 

(fc), was calculated from the Hnear transfer function 
without BAOs based on the fitting formula of Eisenstcin & flu (1998) (see Eq.[29] in their paper). Note that the 
power spectra calculated from the closure theory and RPT are the sum of the leading-order contributions, p(^-'(fc) and 
P'-^^^(fc), not including the higher-order term P^^^^^(fc). For comparison, the one- loop predictions from the standard 
perturbation theory are also depicted as dashed lines. 

On large scales (low-fc), the predictions both from the closure theory and RPT reasonably match the one- loop results 
of standard perturbation theory, as anticipated. This is just the quantitative check for the non-perturbative methods 
discussed in section l4Tl On the other hand, on smaller scales (high-fc), the deviations from the one- loop perturbation 
become manifest and the amplitude of the predictions both from closure theory and RPT is suppressed compared to the 
one-loop predictions. The reduction of the power spectrum amplitude is the natural outcome of the damping behaviors 
appearing in the non-linear propagator and it qualitatively explains the behaviors seen in the N-body simulations (e.g., 
iJeong fc Komatsul 120061 : [Matarrese fc Pietroiiill2007t ICrocce fc Scoccimarrol[2007i ). However, at lower redshifts z ~ 1 
and z = 0.5, the suppression of the amplitude is so significant that the predictions eventually become lower than the 
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Fig. 2. — Power spectrum of density fluctuations Pii(k;z) at z = 1, obtained from the first-order Born approximation to the integral 
solution (see Ea.|64)). The contributions to the total power spectrum are separately plotted as indicated by pW(fc) and pC)(fc) in the 
panel and the total power spectrum, P^^\k) + P^^^\k), is depicted as the dashed lines. Note that in evaluating the power spectrum, the 
approximate solutions for the non-linear propagators G^^^°^ were used. Thick and thin lines indicate the results using the approximate 
solutions G'^^'^°^ from the closure theory and RPT, respectively. 

linear theory prediction. Compared to the prediction from RPT, the suppression of the amphtude is even larger for 
the prediction from the closure theory. 

The strong suppression seen in the low redshifts seems somewhat unphysical, indicating the failure of our present 
analysis. In the panels with redshift z = 1 and z = 0.5 of Figures [3] and |4l the maximum wave number for the 
limitatio n of one-loop perturbatio n theory is estimated according to the condition, A^(fc) = fc'^P(fc)/(27r^) ^0.4, 
given by iJeong fc Komatsul ()2006l ) and is depicted as vertical arrows, above which the perturbation prediction fails 
to reproduce the N-body results within the 1% accuracy. Comparaed to those critical scales, the validity range of 
the non-perturbative predictions seems rather narrow. While this statement seems rather controversial, our present 
analysis actually relies on several approximate treatments. Hence, the inclusion of the higher-order terms such as 
p(m)(fc) might be essential for the correct non-linear behaviors. ICrocce fc Scoccimarrol ( 2007| ) reported that the 
contributions higher than the two-loop correction become important on small scales and the prediction including the 
two-loop correction does improve the prediction, which reasonably agrees with N-body results within an accuracy 
of the N-body calculation. As we noted previously, however, the two-loop correction evaluated by the first-order 
Born approximation is comparable to the higher-order Born approximation coming from the P^^^^\k) term, and at 
a level of the present analysis, it is not clear whether the higher-loop corrections rather than the higher-order Born 
approximation of the one-loop correction are essential or not. Rather, one obvious thing is that a further improvement 
of the approximate treatment is necessary to faithfully predict the non-linear evolution of BAOs. In the line of this, 
self-consistent treatment to solve the evolution equations (|49|) . (fSO]) and ([STjl. would be one plausible approach, which 
we will address this issue in next task. 



6. CONCLUSIONS 

In this paper, motivated by a forthcoming experiment on the precision measurement of the BAOs imprinted on the 
matter power spectrum, a new theoretical tool for predicting the non-linear gravitational evolution of power spectrum 
was presented. In particular, we have applied the non-linear statistical method in turbulence to the cosmological 
perturbation theory and derived a closed set of matter power spectrum. The resultant evolution equations ()49|) - 
(|5T|) with Fourier kernels ((52|) and ((53|) are the non-linear coupled system of integro-differential equations and these 
consistently recover the one-loop results of standard perturbation theory. Further, the exact integral expressions for the 
solutions of closure equations were obtained (see Eqs.[60] and [61] with mode-coupling function [62l). whose analytica l 
expressions coincide with the renormalized one- loop results of the theory developed bv lCrocce fc Scoccimarrol (j2006af l. 
apart from the corrections coming from the vertex renormalization. 

Based on the Born approximation to the exact integral expressions, we next tried to evaluate the non-linear power 
spectrum analytically. Constructing an approximate solution for non-linear propagator, which smoothly matches 
the two asymptotic solutions valid at low-fc and high-fc, power spectrum of density fluctuations was computed and 
the results are compared with those obtained from the RPT. Due to the non-linear damping behavior of the non- 
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Fig. 3. — Ratio of non-linear power spectrum to smoothed linear spectrum, P(fc)/Pno-wigglc('i^)i given at specific redshifts, 2 = 3, 2, 1 and 
0.5. The error bar represents the N-body results taken from ljeong fc Komatsul 1)20061 ). in which different color indicates the results with 
different box size (see their paper in detail). Here, smoothed linear spectra Pno-wigglo{'i'') were calculated from the linear transfer function 
without baryon acoustic oscillation according to the fitting formula of [Eiscnstcin & Hu (1998) (Eq.[29] of their paper). The non-linear 
power spectra are obtained from the first-order Born approximation to the integral solution (Eg. 1641 ). with approximate solutions of the 
non-linear propagator given by closure theory (thick) and RPT (thin). For comparison, one-loop predictions from the standard perturbation 
theory are plotted in dashed lines. Also, in panels with z = \ and 0.5, maximum wave number for limitation of one-loop perturbation is 
indicated by vertical arrows, according to the criterion, A'^(A:) = A;^P(fc)/(27r^) ^ 0.4 (Jeong & Komatsu 2006). 

linear propagator, the resultant amplitude of the power spectra up to the first-order Born approximation is strongly 
suppressed on small scales and this effect becomes significant as decreasing redshift. As a result, predicted spectra at 
high-fc both from the closure theory and RPT fail to reproduce the N-body trends and the inclusion of the higher- 
order corrections is required. Nevertheless, at the intermediate scales, where the damping behavior of the non-linear 
propagator is rather mild, the closure theory and RPT both predict the deviation from the one-loop results of standard 
perturbation theory, which qualitatively agree with the N-body results. 

Although the analytical treatment presented here is still primitive and the range of the applicability is severely 
restricted by the validity of the approximations, the results indicate that our closure theory is a promising non- 
perturbative approach comparable to the RPT and a more elaborated study will provide an accurate prediction for 
non-linear power spectra going beyond standard perturbation theory. In this res pect, a di r ect num erical treatment of 
the closure equations (|49)) - (|5T|) is our urgent task. As it has been reported bv IValageasI ()2007af) . thanks to the full 
numerical treatment, the evolved result of the power spectrum shows several desirable properties and it qualitatively 
matches the N-body trends even at high-fc. In future publication, these points will be further investigated from a 
quantitative point-of-view, particularly focusing on the non-linear evolution of BAOs. 

Finally, one criticism on the present approach including the currently existing non-perturbative methods is that 
the present methodology heavily relies on the single-stream approximation of the Vlasov equation. In a strongly 
non-linear regime at the high-fc region, the shell-crossing eventually occurs and the fluid description would be broken. 
A preliminary investigation suggests that the break down of the single-stream approximation may arise at fc ~ 1 — 
2/iMpc"^ (|Scoccimarrd 120001 1 Afshordil |2007[ ) . Therefore, for the high -fc region beyond the critical scale, the present 
approach cannot be applied and a more delicate treatment based on the Vlasov equation must be developed. This 
issue would be particularly crucial for accurate theoretical predictions to the cosmic shear statistics. 
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APPENDIX 
ONE-LOOP PERTURBATION 

Here, we show that linear plus one-loop power spectra, Pab^\k) +P^t'^\k) + Pal^\k), obtained from the standard 
perturbation theory indeed satisfy the evolution equations given by (j57|) . To do this, we first write down the evolution 

equations for perturbed quantities $a^(fc; 77) in each order: 



A,b(77)<i>p)(fc;77) = 2 



(27r)3 
(27r)3 



Soik - fci - fca) lapqikiM) $(,'^^1; ^7)$^ (^2; v)- 



(Al) 
(A2) 

(A3) 



The above equations are formally solved with a help of the linear propagator, gab{^ ~ v')- For instance, the solution 

(2) 

for second-order quantity $a can be written as 



$(2)(fc;^)= / dr,' g^biv - V') 



d^ki(fk2 



SD{k-k,-k2) 7fcp,(fci,fc2)$W(fci;'y')*i'Hfc2;?7')- 



(A4) 



(27r)3 ^ ' "'^"/v-ii-^y 

Now, let us consider the time evolution of power spectrum, P*^^^^ (fc; 77). Acting the operator Ylabcdkj]) defined by 
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(|20|) on the ensemble average {^^(^\k; rj)^l^\k' ; rj)) , we obtain 

= ($(2)(fc;^) {A,,(r,)cf>f (fc';^)}) + ({A„,(r;)$(2)(fc;,,)} <i>['\k' ; rj)] 
= J dij gai(r] -V) J ^27r)3 J {2tt)^ - ki - ^2) OD{k - fes - fc4) 

x7/r.(fci,fc2)76p«(fe3,fc4)(*['Hfci;'/)$i''(fc2;?yO'&^''(fc3;77)'i>^'nfc4;r?)) 

+ a < — > 6, k < — > fc' ^ , 
= 2 /" drj' gaiiv - v') f d^kT,£k2 Soik' - ki - k.2) S^ik + ki + k2) 

Jrja J 

X 76p9 {ki, k2hirs {ki, fc2)-R^"^ (^1 ; V, '7')-^^"^ (^2 ; V, V') , 
+ ( a< — >b, ki — >fc' ), (A5) 

where, in the last equality, we have replaced the four-point functions with a product of the two-point functions i?^^^^ 
according to Wick's theorem. The quantity i?^^^^ is the linear cross spectrum defined by 

($ii)(fe; ry)$« (fc'; ^')) = {2^f 5n{k + fe') Ril'\k- rj, rj% (r? > ry')- (A6) 

Integrating equation (|A5p over the Fourier mode k' leads to 

^aUv)P^T\k;v) = 2jyv'9ai{v-v') I ^-^^6n{k + k^ + k2) 

X 7;rs(fci, fe2) 76pg(fcl, fc2) i?(,J.^)(fci; ?7, 77')i?^"Hfc2; ?7, r?') 
-t- (a« — >&, fc< — >-fcy (A7) 

Next consider the evolution equation for power spectrum P^^^) (fc; rj). Repeating the similar calculations given above, 
we have 



= (<i>W(fc;^) {AM(^)$i'^(fc';^)}) + ({Aac(77)<i>i^nfc;^)} <i>i'^(fc';^) 

= 2 / d?;' / ^ fl^^x,^^ / ^ fl^'^Nc,^^ <5_D(fc' - fci - fe2) fo(fc2 - fcs - fc4) 

X .ggK'? - V') 76pg(fci, fc2) 7ir..(fe3, fe4)($i'^(fc; v)'^^^\kv, v)'^i}\k3; v'H'\ki; r,') 
a < — f b, k < — > k' 



4 f dr]' ggi{ri - ry') / d^ki(fk2 Soik' - fci - fea) ^^(fe + fci + fea 

X 76pg(fci, fe2)7ir-s(fc, (fc; ^7, ?7')-Rp"^ (fci; v') 

a < — > 6, A; < — > k' 



Here, in the second equality of the above equation, we have used equations (|Aip and (|A3|) , and substituted the formal 

solution Q 
we obtain 



solution (|A4p to rewrite all the perturbed quantities with the linear-order one, (k, 77). Then, integrating over k' , 



^abcdP^f{k;v)=i f drj'g,i{Tj-r,') J fhf^ Soik + k, + k2) 

X Jbp,{ki,k2)jirs{k, fci) r,, TyOi?'^^^^!; ^, 

-h ( a< — >b, k< — > k y (A8) 

Summing up equations (jA7p and (jASp . and using the fact that Habcdiv) Pcd^\k'jV) = Oj we finally arrive at equation 
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INTEGRAL SOLUTIONS 

In this appendix, we show that the integral expressions given in Section r4.2l are compatible with the closure equations 
49|) and (jSOp . Here, we shall particularly focus on the integral expression (|61[) and explicitly derive the closure equation 
50|) from equation ((6T|) . As for the integral expression (l60|) . it is straightforward to show the compatibility between 
equations (^5]) and ([M)) . just repeating the same procedure as presented below. 
Let us consider equation (|6ip and separate the right-hand side of this equation into two terms: 

Rbcif^'-: v') = Gbd{k\Vi m) GceiMv', Vq) Pde{k; %), 

Rbc\k;'r]iV')^ rf»y2 G6d(fch,^?l)Gce(fc|?7',^72)$de(fc;?72,f?l)■ 



'^o 



Vo 



Acting the operator Aabir]) on the above equations, with a help of equation (|5ip. we have 



X Gqi{\k - q\\v, m) Rpr{q;v, Vi) Gsb{k\vi, Vo) Gcd{k\r]' ,T]o) Pbdik-jilo), 



^abiv) Rbc\k;r],r]') = / dr/2 Gce{k\r]' ,r]2) <^ae{k;T]2,r]) 
J no 

+ 4y drji J dri2 J drjs J j^^;^ -/apqiq,k - q)-/irs{-q,k), 

X Gce{k\r]', 772) ^deik; m,Vi) Gqi{\k - q\ \t], 773) Rpr{q; V, V3) Gsd{k\T]3,Vi)- 
Summing up the above two equations, we obtain 



'^abiv) {4c 



(11) I ^ 



dri2 Gceik\T]',ri2)^aeik;ri2,ri) 



Vo 



d^q 
(27r)3 



lapq{q-,k- q)-/irs{-q,k) 



where the bracket [ • • • ] in the second term of right-hand side is rewritten as 



dr]iGqi{\k - q\\r], rji) Rpr{q; V, Vi) Gs&(/clr/i, 770) Gcd{k\ri' , rjo) Pbd{k] m) 



+ dr]i dr]2 / dr]3 Gce{k\r]' ,r]2) ^de{k;ri2,rii) 

Jrio Jrjo J rii 

X Gqi{\k - q\\v,m) Rpr{q;v,V3) Gsd{k\f]3,m) 



drii Gqi{\k ~ q\\r],rii) Rpriq;!^,!!!) 
Gsbik\rii,ria) Gcd{k\r]', 770) Pbdik; 770) 



+ / dr]2 ^773 6(771 -772) G'sd(fcl?7i, 772) Gce(fcl?7', 773) $de(fc;?73,»72) 



drii Gqi{\k - q\\r],rii) Rpriq;'n,r]i) 



Vo 



Rsc{k- rji,Tj') 6(771 - 77') + i?,,(fc; 77', r;i) 6(7;' - 771) 



(Bl) 



Note that, in the second equality, integration variables (771,772,773) has been periodically replaced with (7^2, 7^3, 771) and 
the domain of integral for iji were expanded by introducing the Heaviside step function. On the other hand, with a 
help of the expression ^ab (see Eq.[S5]), the first term in equation (|B1[) becomes 

r' f (Pq 

J ^772 Gce(fc|77', 772) $ae(fc; 7/2,7?) = 2y -J^;^ lapqiq , k - q) -firs{q , k - q) 

r' 

X / dr)iGci{k\r]' ,T]i) Rpr{q;r],rii) Rqs{\k ~ q\;ri,r]i). 
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Now, collecting these results, equation (jBl[) becomes 



no 



[Rscik; ^iVi - V') + Rcs{k- r^' , m) 6(7?' - r;i)} 

r' 

+ 2 dr]i-firs{q,k-q)Gci{k\r]',rii)Rpriq;r],r]i)Rqs{\k-q\;7j,r]i) 

= J lapq{q,k- q) Kcpqi-k,q,k - q;i],r]'), 

with the use of the definition (|53p . This is exactly the same form as in the closure equation (|50p and in this sense, the 
the integral expression (j6f p is compatible with closure equation. 

NON-LINEAR PROPAGATOR 

Here, we summarize the explicit expressions for time-dependent functions a, (3g^d, lg,d and 5 which appear in the 
one- loop propagator G^^' (see also lCrocce fc Scoccimarroll2006bD : 

5 5 



e''/3,(r?)=e-3''/2/3,(77) 

5(^) = |e^''/2_Ie'' + l. 
5 5 

Further, we list the explicit expressions for the scale-dependent functions /, g, h and i (jCrocce fc Scoccimarro|[2"006b() : 



1 



d3q 



504 J (27r)3 fc3q' 



4(fc^-9^)'(2fc2 + 7g2)ln 



k — q 



k + q 



}Po{q), 



1 



d^q 



168 J (27r)3 fcSqS 



k — q 



k + q 



}Po{q), 



d^q 



24 J (27r)3fc3g^ 



d3q 



72 y (27r)3A:3gJ 



+ £(fc2_^2)2(2fc4 + 5fcV + 3g^) In 



k ^ q 



k + q 



k — q 



k + q 



}Po{q), 



'}Poiq). 



Note again that the quantity Poiq) is the linearly extrapolated power spectrum given at the present time. 
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